Many practical systems in physics, chemistry, mechanics, and engineering may involve two timescales,^[1–3] which often behave as a combination of relatively large amplitude and small amplitude oscillations,^[4,5] conventionally denoted by N^K with N and K corresponding to the large amplitude and small amplitude oscillations, respectively. Generally, the system is said to be in quiescent state (QS) when all variables are at rest or exhibit small amplitude oscillations. On the contrary, when the variables exhibit large amplitude oscillations, the system is in spiking state (SP).^[6] Bursting phenomena can be observed when the variables alternate between QSs and SPs in turn. Though a lot of bursting phenomena related to two timescales, such as the relaxation oscillations, have been reported,^[7,8] the dynamics with two timescales received much attention only after the work of Hodgkin and Huxley, who established a simple slow–fast neural model, the oscillations of which agree with the neuron activities.^[9] Divided by the presentation of the slow–fast analysis method by Rinzel,^[10] the related work can be divided into two stages. In the first stage, most of the research focused on the solutions of bursting oscillations, and many asymptotic methods such as quasi-state method and singular perturbation method were proposed,^[11,12] while in the second stage, researchers were devoted to investigate the mechanism of the bursting oscillations by employing the slow–fast analysis method.^[13,14] Two important bifurcations can be observed, occurring at the alternations between QSs and SPs, which can be used to classify the bursting oscillations into different patterns.^[15]

Up to now, most of the results are obtained in the low-dimensional autonomous system, in which QS and SP each has only one form in the bursting oscillations with only co-dimension one bifurcations at the transitions between QSs and SPs. How the high-dimensional dynamical systems with two scales behave remains an open question. Furthermore, when the non-autonomous terms such as periodic excitation involve the vector fields, how to analyze the bursting oscillations needs to be investigated, since the traditional slow–fast analysis method cannot be directly used to explore the bifurcation mechanism.

In this paper, by introducing a periodically changed electrical power source in a coupled BVP oscillator, a five-dimensional dynamical system with periodic excitation is established, which is used as an example to investigate the bursting oscillations in a high-dimensional non-autonomous vector field. By taking suitable parameters with an order gap between the exciting frequency and the natural frequency, different types of bursting oscillations are presented, which may involve multiple forms of QSs and SPs. Furthermore, by regarding the whole exciting term as a slow-varying parameter, the bifurcation mechanism of the bursting oscillations is obtained, which reveals that not only the states of QSs and SPs, but also the bifurcations at the alternations may influence the structures of the bursting attractors.

By introducing a periodically changed electrical power source in the fifth order autonomous system of Wu^[16] which contains two identical oscillators, a circuit model of coupled BVP oscillators with a periodic excitation is obtained, as shown in Fig. 1. A five-dimensional non-autonomous differential equation describing the dynamics of the circuit can be written as

Complicated bursting phenomena can be observed when the exciting frequency is far less than the natural frequency of the system. In this work, we focus on different types of bursting as well as the related mechanisms that would occur when a slow-varying periodic excitation is applied to the coupled BVP oscillators.

Fig. 1.

	Figure Option ViewDownloadNew Window
	Fig. 1. Circuit of coupled BVP oscillators with a periodic excitation.

In order to discover the properties and mechanism of bursting, the stability of the equilibrium points and the related bifurcation are first considered. Here we define the parameter ω ≪ 1 so that an order gap exists between the exciting frequency and the natural frequency, then we investigate the dynamical behaviors of the system with two timescales. Since ω/Ω_N ≪ 1, the trajectories of the state variables may oscillate mainly according to the natural frequency Ω_N, which may be moderated by the small amplitude oscillations with exciting frequency ω. In an arbitrary period T_N starting at τ = τ₀ associated with Ω_N (implying τ ∈ [τ₀,τ₀ + T_N] with T_N = 2π/Ω_N), the whole exciting term w may vary between W_A = A sin(ω τ₀) and W_B = A sin (ω τ₀ + 2π ω/Ω_N), which means that w keeps almost constant during any period T_N corresponding to Ω_N. Therefore, by regarding the whole exciting term w = A sin (ω τ) as a generalized parameter, system (2) could be treated as a generalized autonomous system.

E(x₁(0),y₁(0),x₂(0),y₂(0),x₃(0)) is defined as an equilibrium point of the generalized autonomous system. The characteristic equation of the equilibrium point can be express as

With the existence condition of E and the characteristic equations (3) and (4), the necessary condition of fold bifurcation can be expressed as

Here we only consider two typical cases

Case A By taking k₇ and w as the bifurcation parameters, the bifurcation sets are presented in Fig. 2(a), where the blue and red curves correspond to fold and Hopf bifurcation sets, respectively, which divide the parameter space into several regions associated with different types of dynamical behaviors.

Fig. 2.

	Figure Option ViewDownloadNew Window
	Fig. 2. Bifurcations in Case A. (a) Bifurcation sets on the plane of w – k7. (b) Equilibrium branches as well as the bifurcations for k7 = 1.5. (c) Equilibrium branches as well as the bifurcations for k7 = 4.0.

The equilibrium branches as well as the bifurcations for k₇= 1.5 and k₇= 4.0 are respectively plotted in Figs. 2(b) and 2(c) as examples, in which the solid and dashed curves correspond to stable and unstable equilibrium points, respectively. With the variation of w, the number of equilibrium points may vary between one and three for k₇= 1.5, while it may vary between one and nine for k₇= 4.0. Furthermore, since more equilibrium points are involved in Fig. 2(c), more fold bifurcation points, denoted by A_i (i = 1,2,…,8), can be observed, while only two fold bifurcation points, denoted by LP_i (i = 1,2), exist in Fig. 2(b).

Case B The bifurcation sets for Case B are similar to those for Case A, which we omit here for simplicity. However, great difference between the stabilities of the equilibrium points can be observed. The equilibrium branches as well as the related bifurcations for k₇ = 12.0 are presented in Fig. 3.

Fig. 3.

	Figure Option ViewDownloadNew Window
	Fig. 3. Equilibrium branches as well as the bifurcations for k7= 12.0.

Unlike the situations in Fig. 2(c), though there are the same number of fold bifurcation points, the stabilities of the equilibrium branches are different. One may find that not only the two equilibrium branches E_1± but also E_3± are stable in Fig. 3, which implies that more jumping phenomena can be observed for Case B, leading to multiple stable equilibrium branches to involve in the slow–fast dynamical behaviors.

Here we fix the parameter ω = 0.02 so that an order gap exists between the exciting frequency and the natural frequency and investigate the dynamical behaviors of the system with two timescales. We turn to the two typical cases described above to explore the dynamical evolution of the system.

4.1. Bursting oscillations in Case A

From the bifurcation analysis above, two types of bifurcations, i.e, fold bifurcation and Hopf bifurcation exist in Case A. Furthermore, different numbers of equilibrium points and limit cycles can be observed with the variation of w. Now we fix A = 1.0 and regard k₇ as the bifurcation parameter to explore different bursting oscillations and the related mechanism.

4.1.1. Symmetric subHopf/subHopf bursting

The phase portrait as well as the time history for k₇ = 1.50 is plotted in Fig. 4, from which one can find bursting oscillations with symmetric structure.

Though the period of the movement in Fig. 4 is equal to that of the exciting term at T = 2π/ω ≈ 314.16, in the neighborhoods of the two equilibrium points E₁₊ and E₁₋, relatively large-amplitude oscillations on the trajectory can be observed, with the frequency of approximately 2.50. Accordingly, the trajectory can be divided into four segments, corresponding to two quiescent states and two repetitive spiking states.

Fig. 4.

	Figure Option ViewDownloadNew Window
	Fig. 4. Bursting oscillations for k = 1.50: (a) phase portrait, (b) time history, (c) local enlarged time history of (b).

Fig. 5.

	Figure Option ViewDownloadNew Window
	Fig. 5. Transformed phase portrait on the plane of (w,x1) for k7 = 1.5.

In order to reveal the mechanism of the bursting oscillations in Fig. 4, we present the overlap of the transformed phase portrait with equilibrium branches in Fig. 5. According to the time history shown in Fig. 4(b), it is easy to find that the motion direction of the trajectory is clockwise. If the trajectory starts at the point B in Fig. 5, it moves almost strictly along the equilibrium branch E₁₊, which is in QS. When the trajectory passes through point subH₁, E₁₊ loses stability via subHopf bifurcation. It should be pointed out that the trajectory does not jump to the stable equilibrium point E₁₋ immediately with the attraction of E₁₋, instead it moves along the unstable equilibrium branch E₁₊ for a while because of inertia. This phenomenon is called the delay behavior.^[17] Then the trajectory jumps to E₁₋ and oscillates around it, leading the system from QS to SP. With slow variable w continuously increasing, the trajectory gradually converges to the stable equilibrium point E₁₋ and transforms to QS again. Until w reaches to its maximum 1.0, the movement turns to the left. Due to the symmetry of the vector field, the left part of the trajectory is similar to the right part. Till the trajectory returns to point B, a periodic bursting is completed. Obviously, the repetitive spiking states are created since slow variable w passes through the subHopf bifurcation point. So such a bursting can be referred to as symmetric subHopf/subHopf bursting of point–point type.

4.1.2. Symmetric subHopf/fold-cycle bursting

With the increase of k₇, different types of bursting attractors can be obtained, an example is shown in Fig. 6 with k₇ = 4.0.

Different from the situation in Fig. 4, though the trajectory in Fig. 6 can also be divided into four segments, each repetitive spiking state is composed of two parts associated with relatively large-amplitude oscillations SP₁ and relatively small-amplitude oscillations SP₂ (see Fig. 6(c)). Note that the period of the bursting oscillations in Fig. 6 is the same as that in Fig. 4. However, when the trajectory arrives at the transition points from QS to SP, for example, at the point S in Fig. 6(c), instead of jumping from the equilibrium branch E₁₊ to E₁₋, it oscillates almost along the stable limit cycle of the generalized autonomous system CY shown in Fig. 6(d), resulting in large-amplitude oscillations SP₁. After twice around CY, the trajectory settles down to the equilibrium branch E₁₋, leading to small-amplitude oscillations SP₂.

Fig. 6.

	Figure Option ViewDownloadNew Window
	Fig. 6. Bursting oscillations for k7 = 4.0: (a) phase portrait, (b) time history, (c) local enlarged time history of (b), (d) stable limit cycle in the generalized autonomous system.

The equilibrium branches and the limit cycles of the generalized autonomous system with the variation of w are shown in Fig. 7(a), in which the red and blue curves respectively correspond to stable and unstable limit cycles, while LPC₊ and LPC₋ represent the fold bifurcation points of the limit cycle. To clarify it, we present the stable equilibrium points and the limit cycles with different w in Table 1. From the transformed phase portrait in Fig. 7(b), the associated bifurcation mechanism of bursting can be revealed explicitly, in which the solid and hollow cycles correspond to stable and unstable limit cycles, respectively. A trajectory starts at point B, it moves along the equilibrium branch E₁₊, which is in QS. When it passes through subH₁, E₁₊ loses its stability. The trajectory also generates the delay behavior, moving along the unstable equilibrium branch E₁₊ until w = 0.118. It must be especially emphasized that the stable limit cycle and equilibrium branch E₁₋ coexist in the same phase space, but the trajectory is in the attraction domain of the limit cycle, and therefore, the trajectory appears in large-amplitude oscillation and transforms to SP₁ from QS. Until the stable limit cycle disappears at LPC₋, the trajectory oscillates around the stable equilibrium branch E₁₋, entering into relatively small-amplitude oscillations SP₂. With the increase of w, the trajectory gradually converges to E₁₋, leading to QS from SP₂. When the variable w reaches the maximum value 1.0, the trajectory moves to the left. Due to the symmetry of the vector field, similar dynamical behaviors appear with the variety of w until it reaches to point B, completing a periodic bursting.

Fig. 7.

	Figure Option ViewDownloadNew Window
	Fig. 7. Bifurcation mechanism of bursting for k7 = 4.0: (a) equilibrium branches and limit cycles, (b) transformed phase portrait on the plane of (w,x1).

Table 1.

Table 1.

Table 1. Equilibrium points and limit cycles with variation of w for k7 = 4.0. . w Stable equilibrium points Limit cycles stable unstable (−∞, −0.4619] E1+ No No [−0.4619, − 0.1386] E1+ Yes Yes [−0.1386,0.1386] None Yes No [0.1386,0.4619] E1− Yes Yes [0.4619, +∞) E1− No No

Table 1. Equilibrium points and limit cycles with variation of w for k7 = 4.0. .

In this case, SPs and QSs are connected by the subHopf bifurcation of the equilibrium points and the fold bifurcation of the limit cycle. So, we refer to such a bursting attractor as symmetric subHopf/fold-cycle bursting.

4.2. Bursting oscillations in Case B

The attractors of bursting oscillations in Case B are different from those in Case A, since the properties of the equilibrium points are different. Figure 8 plots details of the bursting oscillations for k₇ = 12.0, A = 5.0 as an example, in which the trajectory can be divided into eight segments, corresponding to four quiescent states and four repetitive spiking states, denoted by SP_i and QS_i (i = 1,2,…,4) in Fig. 9.

Fig. 8.

	Figure Option ViewDownloadNew Window
	Fig. 8. Bifurcation mechanism of bursting for k7 = 12.0: (a) portrait phase of three-dimensional, (b) portrait phase of two-dimensional, (c) time history, (d) local enlarged time history of (c).

Since w is a slow variable, combining with the equilibrium branches shown in Fig. 3, the bifurcation mechanism of bursting can be obtained upon the transformed phase plotted in Fig. 9. A trajectory starts at point B, it moves almost along the equilibrium branch E₁₊. The trajectory is in QS₁ until it arrives at point A₁, where a fold bifurcation takes place. Although two attractors E₃₊ and E₁₋ coexist in the same phase space, the trajectory is in the attraction basin of E₃₊, and thus the trajectory jumps to the stable equilibrium point E₃₊, leading the system to SP₁. Then the trajectory oscillates and converges to E₃₊ with the increase of w, which causes SP₁ to settle down to QS₂. When the trajectory arrives at another fold bifurcation point A₅, E₃₊ disappears with fold bifurcation, and the trajectory jumps to the stable equilibrium point E₁₋ to form SP₂. Then the trajectory gradually converges to E₁₋ and transforms to QS₃. With further increase of w, it reaches the maximum value 5.0, and the trajectory moves in the left direction under the effect of the external excitation. Because of the symmetry of the vector field, when the trajectory returns to point B, one period of the bursting oscillations is finished which contains four transitions between QSs and SPs.

In this bursting oscillation, the trajectory undergoes four QSs and four SPs during each period, and the transition between QSs and SPs is caused by fold bifurcations. Therefore we can refer to such a bursting attractor as symmetric double-fold/fold bursting of point–point type.

Fig. 9.

	Figure Option ViewDownloadNew Window
	Fig. 9. Transformed phase diagram on the plane of (w,x1) for k7 = 12.0.

Different types of periodic bursting oscillations for the coupled BVP oscillators with external excitation can be observed when an order gap exists between the external exciting frequency and the natural frequency. By introducing the concept of the generalized autonomous system in which the whole exciting term is regarded as a slow-varying parameter, the bifurcation forms as well as the critical conditions of the equilibrium points are presented. Based on the bifurcation analysis, two cases are chosen to show the typical bursting oscillations of the system. Under several parameter conditions corresponding to different bifurcation forms involving the oscillations, different forms of bursting oscillations such as symmetric subHopf/subHopf bursting, symmetric subHopf/fold-cycle bursting, and symmetric double-fold/fold bursting can be observed, the mechanism of which is explored by employing the transformed phase portraits. It must be pointed out that since multiple attractions may coexist in a high-dimensional system, the approaching of the trajectory depends on the attractors that the trajectory settles down to, which may determine the types of QSs and SPs. Furthermore, the delay behavior can be observed at the subHopf bifurcation point which connects the QS and SP, while other bifurcation points agree well with the turning points on the trajectories.